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If the Pioneer anomaly has a gravitational origin, it would, according to the equivalence principle, 
distort the motions of the planets in the Solar System. Since no anomalous motion of the planets has 
been detected, it is generally believed that the Pioneer anomaly can not originate from a gravitational 
source in the Solar System. However, this conclusion becomes less obvious when considering models 
that either imply modifications to gravity at long range or gravitational sources localized to the 
outer Solar System, given the uncertainty in the orbital parameters of the outer planets. Following 
the general assumption that the Pioneer spacecraft move geodesically in a spherically symmetric 
spacetime metric, we derive the metric disturbance that is needed in order to account for the Pioneer 
anomaly. We then analyze the residual effects on the astronomical observables of the three outer 
planets that would arise from this metric disturbance, given an arbitrary metric theory of gravity. 
Providing a method for comparing the computed residuals with actual residuals, our results imply 
that the presence of a perturbation to the gravitational field necessary to induce the Pioneer anomaly 
is in conflict with available data for the planets Uranus and Pluto, but not for Neptune. We therefore 
conclude that the motion of the Pioneer spacecraft must be non-geodesic. Since our results are model 
independent within the class of metric theories of gravity, they can be applied to rule out any model 
of the Pioneer anomaly that implies that the Pioneer spacecraft move geodesically in a perturbed 
spacetime metric, regardless of the origin of this metric disturbance. 

PACS numbers: 04.80. Cc, 95.10.eg, 95.30. Sf, 95.55.Pe, 96.30.-t 



I. INTRODUCTION 

A. Brief summary of the Pioneer anomaly 

Precision tracking of the Pioneer 10 and 11 spacecraft 
over more than 30 years has revealed an anomalous fre- 
quency drift of the Doppler radar signals transmitted by 
the spacecraft [3, 0, 0, 0, 0]. After careful analysis 
of potential sources of this frequency drift, it has been 
concluded that it most likely is caused by an anoma- 
lous acceleration of both spacecraft, directed towards the 
Sun. Moreover, the anomalous acceleration seems to stay 
constant over long distances, and it has been estimated 
to have very similar values for both spacecraft, about 
8 x 10 _10 m/s 2 over distances between 20 and 70 AU from 
the Sun Q. The data presented by Anderson et. al. @,IH 
show that the anomalous acceleration became noticable 
beyond the orbit of Saturn, i.e. around 10 AU. It quickly 
grew to a near constant value of around 8 x 10~ 10 to/s 2 
directed towards the Sun, and has remained more or less 
constant beyond 20 AU. As the Pioneer anomaly remains 
unaccounted for, the question has been raised whether 
the anomaly could originate from unknown physics. In- 
deed, a number of proposals for possible physical expla- 
nations have been made. These range from conventional 
physics, such as anisotropic heat dissipation 0,01; drag 
from interplanetary dust |8fl or gravitational attraction 
from Kuiper-belt objects [3], to new physics, such as 
scalar fields [10], brane- world scenarios [ll[, dark matter, 
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Modified Newtonian Dynamics (MOND) [H, E3], mod- 
ifications to the Newtonian gravitational potential [uj |. 
modifications to General Relativity [l5| , and other alter- 
native gravitational theories, see refs. 0, 0] and refer- 
ences therein. 

A natural question to ask is whether the effect could 
have a gravitational originj3(3], either originating from an 
unknown gravitational source present within the Solar 
System, or being a manifestation of long-range modifi- 
cations to gravity, i.e. arising from deviations in the 
vacuum gravitational field from the field predicted by 
General Relativity, where these deviations become ob- 
servable only over longer distances. If the latter were 
found to be true, it would imply the need for a revision of 
established gravitational theory, a possibility that, how- 
ever remote it might be, by itself explains the attention 
the Pioneer anomaly has received within the fundamen- 
tal physics community. If the Pioneer anomaly were in- 
deed a gravitational effect, one would, according to the 
equivalence principle, expect such a disturbance in the 
gravitational field to also affect the motions of the plan- 
ets. The question of whether the Pioneer anomaly could 
be attributed to a gravitational source in the Solar Sys- 
tem was considered by Anderson et. al. The orbits 
of the inner planets are determined with very high pre- 
cision, and no anomalous motion of these planets have 
been observed, so if a gravitational source were to ex- 
plain the Pioneer anomaly, it could not be present in 
the inner Solar System. However, because the Pioneer 
effect first became visible as the spacecraft entered the 
outer Solar System, their analysis does not exclude the 
possibility of a gravitational source localized to the outer 
regions of the Solar System, beyond the orbit of Saturn. 
In this case, one would expect the disturbance in the 
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gravitational field created by such a source to only affect 
the motions of the three outer planets, which orbit the 
Sun at distances between 20 and 40 AU, but it would 
have no observable effect on the motions of the inner 
planets. Since the orbits of the three outer planets have 
been determined by optical methods alone [H, [TtJ , the 
uncertainty in their orbital parameters is much higher 
than for the other planets. Recently, Page, Dixon and 
Wallin Izzo and Rathke [llj and Iorio and Giudice 
[20| examined the question whether the Pioneer anomaly, 
given that it is a gravitational effect, could have an ob- 
servable effect on the orbits of the outer planets. Page 
et. al. concluded that, if the Pioneer anomaly has gravi- 
tational origin, it would not have an observable effect on 
the motions of the outer planets, given the present accu- 
racy of planetary data for these planets. Their conclu- 
sion was based on the fact that only optical observations 
have been made of Uranus, Neptune and Pluto [lfl E3] 
and that the a priori uncertainty of the optical measure- 
ments for the three outer planets is substantial, resulting 
in rather large uncertainties in their orbital parameters. 
Such uncertainties could possibly mask the disturbances 
to the planetary orbits caused by an anomalous gravi- 
tational field. However, no calculations were shown to 
substantiate their conclusion. Izzo and Rathke, on the 
other hand, came to a different conclusion. They con- 
cluded that the Pioneer anomaly could not have gravi- 
tational origin, because it would be in conflict with the 
ephemerides of Uranus and Neptune. They used para- 
metric constraints set on a parametrization of deviations 
from Newtonian gravity in terms of an effective reduced 
Solar mass Ho(r) [2l|. Applying the constraints to the 
planets Uranus and Neptune, they concluded that the de- 
viations to the Newtonian gravitational potential implied 
by the Pioneer anomaly would lead to changes in the ef- 
fective reduced Solar mass that exceeded the limits set 
by Solar System measurements by 1-2 orders of magni- 
tude. Iorio and Giudice [20] came to the same conclusion, 
providing an assessment of the perturbations to the or- 
bits of the outer planets based on the Gauss equations. 
The Gauss equations estimate the effect of a perturba- 
tion in the gravitational field in terms of the time rate of 
change it causes on the osculating orbital elements. Iorio 
and Giudice found large secular effects on the argument 
of the perihelion and mean anomaly for all three outer 
planets. 



The possibility that the Pioneer anomaly might indi- 
cate a deviation in the vacuum gravitational potential 
from the one predicted by General Relativity has also 
been studied. A possible Yukawa modification to the 
Newtonian gravitational potential has been dismissed be- 
cause of constraints from planetary data @,[H[. However, 
modified gravitational theories have been proposed that 
are claimed to account for the Pioneer anomaly without 
viola ting parametric constraints inferred from planetary 
data [3]. 



B. What is the appropriate question? 

As we can see from this brief review of the subject 
(which by no means is complete, for a more complete re- 
view see refs. 0, Q), the number of models that have 
been proposed to explain the Pioneer anomaly is already 
significant, and growing fast. Fitting a model to paramet- 
ric constraints is useful for dismissal and initial scrutiny 
of models, but a closer examination of the models that 
pass the parametric constraints is in place. When look- 
ing at the models aiming to explain the Pioneer anomaly, 
they are different, but many of them are based on the as- 
sumption that the motion of the Pioneer spacecraft is 
geodesic. The question whether the Pioneer effect could 
have gravitational origin, either being caused by an un- 
known gravitational source localized in the outer regions 
of the Solar System or arising from long-range devia- 
tions in the gravitational field from the vacuum field pre- 
dicted by General Relativity, can be reformulated as a 
question of whether the Pioneer spacecraft move geodesi- 
cally. Within the class of metric theories of gravity, this 
question is model independent. As will be seen, it can 
be satisfactorily answered by giving quantitative, model- 
independent estimates derived from the basic equations 
of motions, of the perturbations to the orbits of the outer 
planets that would result from such an effect. These es- 
timates could then be compared with observations, thus 
providing a more thorough test of models against obser- 
vations. To settle this question would be of significance 
for the substantial model building effort targeting the Pi- 
oneer anomaly. Furthermore, if it were possible to prove 
that the Pioneer spacecraft moved non-geodesically, it 
would imply that we could rule out a large class of mod- 
els, including several of the models that have been pro- 
posed to explain the Pioneer anomaly. 



C. And how can we answer it? 

The purpose of this paper is to derive an answer to 
the question posed above; whether the Pioneers move 
geodesically. We will perform a quantitative, model- 
independent assessment of the perturbations to the mo- 
tions of the outer planets that would arise from a dis- 
turbance to the spacetime metric of the kind required to 
account for the Pioneer anomaly. The rationale for ana- 
lyzing the motions of the outer planets is, in addition to 
the obvious reason that they orbit the Sun at distances 
at which the Pioneer anomaly has been observed, the 
relatively large uncertainty in their orbital parameters 
El- Therefore, the likelihood of our result will to a 
large extent depend on this uncertainty. 

The procedures that are currently applied for calcu- 
lating the ephemerides of the planets use available the- 
ory, taking into account known gravitational sources, and 
carefully fitting the computed orbits to available obser- 
vational data [H, [TtJ ■ One would therefore expect that 
any disturbance to the gravitational field within the solar 
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system would, if large enough, over time generate resid- 
uals in the position measurements of the planets as they 
would tend to drift from their expected positions. The 
question then is, what residuals would arise if we were to 
explain the Pioneer anomaly as the result of a perturba- 
tion to the gravitational field in the outer solar system, 
and how do they compare with the observed residuals? 

Whether or not this metric perturbation is induced 
by a real gravitational source or is just a manifestation 
of a long-range deviation in the vacuum gravitational 
field from the field predicted by General Relativity is ir- 
relevent for our study. In order to determine whether a 
perturbation to the gravitational field has an observable 
effect on planetary orbits, we will have to simulate the 
effect that the gravitational perturbation would have on 
the planetary orbits (what we will refer to as the gross 
perturbations or gross effects), then compute the simu- 
lated ephemeris as the solution to the unperturbed equa- 
tions of motion that provides the best fit to the perturbed 
solution, compute the simulated residual as the difference 
between those two solutions and finally compare the sim- 
ulated residuals with the real residuals. 

An assessment similar to our study was recently done 
by Iorio and Giudice [2(J, based on the Gauss equations 
for orbital motion. The Gauss equations estimate rates 
of change of orbital elements, an approach that makes 
the effects easy to understand. On the other hand, the 
orbital elements are not directly observable, and the esti- 
mation of observables from them is not straight forward. 
In their work, they investigated the gross orbit pertur- 
bations caused by a constant, anomalous acceleration of 
the Pioneer kind. They found large secular effects on 
the argument of the perihelion and mean anomaly for all 
three outer planets. As will be shown later, such large 
secular effects will largely be canceled by the ephemeris, 
because secular effects are lowest-order effects that to a 
large extent will be matched by a best fit solution to the 
unperturbed equations of motion. 

Our approach is to work from the basic equations of 
motion, thus working with the real orbit of the body, not 
with approximate orbital elements. Relating our results 
to measurements is only a matter of making coordinate 
transformations between the reference frame of the calcu- 
lation and the reference frame of the observer. We apply 
the hypothesis that the Pioneer spacecraft move geodesi- 
cally through space, i.e. that the only forces acting on 
them are gravitational and that the anomalous accelera- 
tion is the effect of disturbances to the spacetime metric. 
Furthermore, we assume a spherically symmetric gravi- 
tational field that can be described within the framework 
of a metric theory of gravitation. The spacecraft move 
approximately in opposite directions and their trajecto- 
ries have both minor inclinations relative to the ecliptic 
plane [3], so a spherically symmetric gravitational field 
at large distances from the Sun should be a reasonable 
assumption. Since we are assuming a metric theory of 
gravitation, the equivalence principle will be valid within 
the theory [H, H3|- At the outset, we assume nothing 



about how the metric couples to gravitational sources or 
of the nature of any gravitational source, so our results 
will be model independent in the sense that they will 
only depend on the assumptions of geodesic motion in a 
static spacetime metric and spherical symmetry. We will 
work in the weak field limit, and in this limit, a metric 
theory of gravitation can be parametrized using the PPN 
formalism [22, [23j. Thus, our results will be valid for any 
PPN-parametrizable gravitational theory. 

We are only interested in dominant effects and will 
work perturbatively to lowest order in the non-relativistic 
weak field limit. To lowest order, independent pertur- 
bations decouple, so we may disregard other effects that 
perturb the planetary orbits, like disturbances from other 
planets or satellites. This simplification allows us to treat 
the problem analytically by doing perturbative expan- 
sions around idealized analytical solutions to the plane- 
tary orbits. 



D. Article overview 

In sectionHIlwe derive and solve the perturbative equa- 
tions of motion for Newtonian gravity. We make a quan- 
titative assessment of the gross perturbations in order to 
get a better understanding of the real effects of the gravi- 
tational perturbation. In this section we also account for 
the accuracy of the perturbative expansions applied. 

Then, in section HTT1 we make the generalization of this 
procedure to any metric theory of gravitation. First, we 
derive the unique form that the perturbation to the grav- 
itational field must have in order to be able to account for 
the Pioneer anomaly. From the equations of motion, it is 
possible to uniquely derive the signature of the gravita- 
tional potential that would cause a constant anomalous 
acceleration. Having uniquely constrained the spacetime 
metric, we are able to derive the perturbative equations 
for geodesic motion and find the general solution to them. 

In section Hv| we cover in detail the analysis method 
that will be employed in order to determine whether the 
perturbations computed in sections [II] and IIIII are ob- 
servable or not. 

Finally, in section [Vl we apply this analysis method to 
the results of sections Ull and HTT1 Here, we first derive the 
simulated ephemeris, which is the solution to the unper- 
turbed equations of motion that best fits the perturbed 
solution. Having obtained the simulated ephemeris, we 
then compute the simulated residuals as the difference 
between the observable values of the perturbed orbit 
and the simulated ephemeris. Finally, we are able to con- 
duct the result analysis described in Section llVl where we 
compare the simulated residuals with the real residuals. 
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II. PERTURBATIONS TO PLANETARY 
ORBITS - NEWTONIAN GRAVITY 

We will start our analysis with the simplest case, which 
is that of a static, spherically symmetric, Newtonian 
gravitational field. 



A. Equations of motion 

For an object moving freely in a static, spherically 
symmetric gravitational field, motion is planar. In or- 
der to describe the motion of the object, let us therefore 
introduce in-plane coordinates (r,cp), where r is radial 
distance from the central mass to the object, and ip mea- 
sures the longitude of the object in the orbital plane, 
relative to a fixed reference direction. In the case of 
planetary motion, we will define ip as the argument of 
latitude, which is positive in the direction of movement 
and measured relative to the ascending node, the point 
where the body ascends through the ecliptic plane. The 
equations of motion can easily be reduced to the follow- 
ing form, expressed in terms of the resiprocal distance 
u=l/r 0: 



rn 



dcp(t) 



dt 



d 2 u(ip) 
dip 2 



i 2 u 2 ' 



(1) 

(2) 



where the constant t is the specific angular momentum 
(angular momentum per unit mass) and !F(u) is the force 
per unit mass. Any solution to these equations have two 
constants of motion: The angular momentum L and the 
total energy E, which is the sum of kinetic and potential 
energy. For a Newtonian gravitational field surrounding 
a central mass M, the force !F(u) has the familiar inverse 
square form: Tn = — where G is the gravitational 
constant. In this case, eq. [2] takes the form 



d 2 u(ip) 
dip 2 



GM 



(3) 



The bound solutions to eq. [3] are the familiar Keple- 
rian orbits with elliptical form: 



ro(<P) = 



1 



a(l-e 2 ) 



uo(ip) 1 + e cos(tp — uj) ' 



(4) 



where a is the semi-major axis of the orbit, e is the 
eccentricity and ui is a constant, the argument of periapsis 
of the orbit. Notice that the angular difference ip — li is 
the true anomaly of the orbit, denoted v. 



B. Pioneer perturbation 



Let us perturb the Newtonian gravitational field with 
an anomalous, constant radial accelerations a p — 8.7 x 



10" 
the form 



I s 2 directed towards the Sun. The force T takes 



GM 



(5) 



Assume that uq is a reference solution to the unper- 
turbed equations of motion, of the form given by eq. |H 
Define ujq to be the argument of periapsis of the reference 
orbit. Then, let us define the perturbation Su — u — Uq, 
where u is a solution to the perturbed equation of motion, 
eq. El with the force T defined by eq. 03 

From eq. we obtain the following equation of motion 
for the perturbation Su: 



d 2 Su 



1 



2k 



Su 



dip 2 \ ' (1 + ecosvf J ' a(l - e 2 )(l + ecosi/) 2 

(6) 

where the constant n = ^-(1 — e 2 ) is proportional 
to the ratio between the anomalous acceleration a p and 
the mean Newtonian acceleration ajv = GM/a 2 , and we 
have introduced v = <p — uio for notational convenience. 
Evaluating k at heliocentric distances relevant for this 
study yields n < 10 -4 . The presence of K on the left hand 
side of eq. [6] gives rise to a small perihelion precession 
of the orbit. However, the contribution from the K-term 
is small compared with the total perturbation. This can 
easily be seen by taking the circular limit of eq. \E[ in 



which case Su 



a(l+2ft) 



For eccentricities e < 1, 



this will still be the case. Subsequently, we will disregard 
the angular dependency of the K-term on the left hand 
side of eq. [6] and use the approximate equation j3]J 



d 2 5u , „ . . k 

— s- + (1 + 2k) Su = j • 

dp a(l — e 2 )(l + ecosf) 



(7) 



Initial conditions for the perturbation are chosen by 
demanding that the solutions u and uq have the same 
angular momentum I (Si — 0). Furthermore, we will 
demand that the perturbations Su and Sip vanish iden- 
tically in the limit K — > 0, i.e. the perturbed orbit ap- 
proaches the reference orbit as the perturbation vanishes. 
It should be noted that there is no loss in generality here; 
since there are no additional degrees of freedom, the free- 
dom in the choice of initial condition for the perturbed 
solution u lies in the freedom of choosing the reference 
orbit uq. If we do a series expansion in eccentricity e 
of the right-hand side of eq. 0, the equation can be inte- 
grated analytically. It turns out that we can get sufficient 
accuracy by expanding to third order in eccentricity, cf. 
Section III CI below. In order to simplify the expressions 
and maintain focus on the dominant terms, we will leave 
out the higher order terms in the remainder of this sec- 
tion. The reader may find the corresponding third order 
expressions in Appendix [Aj It should be pointed out, 
though, that all perturbative results of the present pa- 
per were derived expanding to third order in eccentricity. 
The general solution of eq. [7] is, to lowest order in e and 
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k: 

. K e 

Su = cos v 

a a 

+ Asin(Vl + 2ku) +Bcos{y/l + 2kv) + 0{e 2 ), (8) 

where A and B are constants of integration. Notice 
that v is the true anomaly of the reference orbit. De- 
manding Su = for k = fixes A and B to take the 
values A — and _B = e/a + 0(e 2 ). The perturbation 
then takes the following form to lowest order in e and k: 

Su = ^(l-ej/sinz/) + 0(e 2 ). (9) 

This gives the perturbation in radial distance, Sr: 

Sr = — = — Ka(l - 2ecos^ - et/sinz/) + 0(e 2 ). (10) 

Thus, the gross effect of a Pioneer type perturbation 
to the Newtonian gravitational field is a reduction in the 
orbital radius of the order of na ~ a— . It should be 

ajv 

stated here that the radial perturbation Sr, as given by 
eq. [TOl differs from the expression given by eq. 58 of 
Anderson et. al Q. They provide an expression for Sr 
that does not vary with the longitude ip, and therefore 
can not be a solution to the equation of motion, eq. [7] 
Eq. [10] coincides with eq. 58 of 0] in the circular limit, 
though. Notice that from eqs. [9] and [10] we see that 
Sr ~ nro 2 /a. This implies that the more eccentric orbits 
should show much stronger variation in the perturbation 
over the orbital period than less eccentric orbits. In the 
circular limit, the perturbation in orbital radius is con- 
stant, Sr e= Q = —na. 

Next, let us turn our attention to eq. [TJ We can use 
it to obtain an equation for the perturbation Sip. Using 
that Si = 0, and that ^ = = we get the 

following equation for the angular perturbation Sip: 

^£ = 2^ (11) 

dip ito 

Using eq. we get 

— — =2K(l-ecos^-ei/sini/)+0(e 2 ). (12) 
dip 

We can use eq. [J3] to obtain an expression for the 
perturbation to the angular velocity if: 

dSip dip y 

Sip = — — = 2k (1 — ecos v — evsmv) if + U(e ) 

dip dt 

(13) 

Thus, the angular velocity increases by an amount of 
the order of 2mp, which implies a perturbation to the 
mean motion of the same order. This is consistent with 
the conservation of angular momentum £ = r 2 ip, which 
implies that Sip = —2^-ip. As the planet is drawn closer 
to the Sun due to an increased gravitational pull, the 



angular velocity increases in order for the angular mo- 
mentum to remain constant. 

Integrating eq. [12] and demanding that Sip vanishes for 
k = 0, we obtain the perturbation Sip: 

Sip = 2k(v + evcosv - 2esin^) + 0(e 2 ) (14) 

This is the gross perturbation to the argument of lat- 
itude that should be expected from a Pioneer type per- 
turbation to the Newtonian gravitational field. The per- 
turbation contains two terms: The first term is linear in 
v, and represents a big secular effect, Sip^ — 2kv. It 
comes from an increase in the mean motion balancing 
the overall negative perturbation to the orbital radius. 
The other term is a smaller ripple on the secular term. 

A precise statement should be made regarding the use 
of the variable <p in the above expressions: ip{t) is the 
argument of latitude of the reference orbit at a particular 
time t and therefore labels the perturbations uniquely in 
time: Sr(t) = Sr{ip{t)) and Sip(t) = Stp(ip(t)). 



C. Accuracy of the perturbative approach 

A few comments are in place regarding the accuracy 
of the perturbative approach used in Section III Bt When 
integrating eqs. [7] and [12] in order to compute the per- 
turbations Su and Sip, care must be shown when doing 
series approximations with respect to the eccentricity of 
the right-hand sides of the equations. Pluto has a rel- 
atively high eccentricity of ~ 0.25. Comparison of ex- 
act numerical with pertubative solutions revealed that 
we should expand to 3rd order in eccentricity in order to 
get sufficient accuracy in the perturbations to the orbit 
of Pluto. A comparison of results from numerical inte- 
grations over a complete period of the exact, non-linear 
equation of motion, eq. [21 with the perturbative New- 
tonian equations of motion, eq. [7] shows that the per- 
turbative solutions for the radial perturbation Su have 
acceptable accuracy for all three outer planets when we 
expand the right-hand side of the equation to third order 
in eccentricity. The relative errors A(Sr)/Sr in the radial 
perturbation Sr were 0.02% and 1.9% for Uranus and 
Pluto, respectively. Although these are still big errors in 
absolute terms due to the large distances involved, they 
will not have significant impact on our results, because 
we are seeking dominant effects only. A similar compar- 
ison of perturbative and numerical solutions to eq. [H] 
shows that the relative errors A(Sip)/Sip were 5 x 10~ 6 
and 0.5% for Uranus and Pluto, respectively. Finally, the 
parameter that defines the perturbation to the gravita- 
tional field is the parameter k, which is of the order of 
the ratio between the anomalous acceleration a p and the 
mean Newtonian gravitational acceleration: k ~ a p /a,N- 
k has a small value for all three planets, of the order of 
10 -4 for Pluto, and less for the other two. A comparison 
of exact and perturbative solutions reveals that a first 
order perturbative expansion in k gives good accuracy. 
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D. Assessment of gross effects for the three outer 
planets 

Let us apply the results of the previous section to the 
three outer planets. Eq. [10] gives the perturbation Sr to 
orbital radius as a function of true anomaly v = ip — luq. If 
we evaluate the e° perturbation Sr^ = — na for the three 
outer planets, we get values of —155,000 km, —598,000 
km and —1,359,000 km, respectively, for Uranus, Nep- 
tune and Pluto. Figure Q] plots 5r(v) for the three outer 
planets. As we can see from the figure, Pluto shows a 
much stronger variation in Sr over the orbital period than 
the other two planets. This is due the relative high eccen- 
tricity of Pluto's orbit (~ 0.25) and confirms the strong 
coupling of Sr to the orbital radius mentioned in the pre- 
vious section. 



Perturbation to orbital radius 
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FIG. 1: Gross perturbation to orbital radius (10 3 km) over a 
complete period for the three outer planets 
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Perturbation to angular velocity 
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FIG. 2: Gross perturbation to angular velocity in arc sec/year 
over a complete period for the three outer planets 
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The perturbation Sip to the angular velocity is given by 
eq. and the e° perturbation is Sip^ = 2nip. Evaluat- 
ing it for the three outer planets gives values of 1.7,2.1 
and 2.4 arc sec/year for Uranus, Neptune and Pluto, re- 
spectively. <fi{v) is plotted in Figured The variation in 
the perturbation Sip over the orbital period is shown in 
Figure[2l As we can see, Pluto has the strongest variation 
in angular velocity, which is because of the proportional- 
ity Sip oc —5r/r 3 . Sip increases as \fa with orbital radius, 
and we can clearly see from Figure[2]that there is stronger 
variation in the angular velocity perturbation for orbits 
with higher eccentricity. 

The perturbation to the argument of latitude, Sip, is 
given by eq. HH The e° term of the perturbation to 
the angular velocity Sip implies a large secular increase 
of Sip, which is apparent in the plot of Sip in Figure 03 
For Uranus and Neptune, the secular effect is completely 
dominant, and Sip increases almost linearly with ip. For 
Pluto, however, the big variation in Sip that can be seen 
in Figure [2] makes Sip for Pluto deviate significantly from 
a linear increase. 

If these perturbations were real, one would expect that 
by fitting a Keplerian orbit to them that has a slightly 



FIG. 3: Gross perturbation to the argument of latitude ip (arc 
sec) over a complete period for the three outer planets 



higher mean motion, one would be able to cancel the 
secular increase. This can easily be seen by the following 
example: If the mean motion was increased by an amount 
equal to the secular term Sip^ = 2np, we would expect 
the resulting residuals to show no secular increase. Figure 
|4] shows the net perturbation Sip — Sip^ for Uranus and 
Neptune after removing the secular term Sip^ . As would 
be expected, the net perturbations are much smaller, of 
the order of a few arc seconds. Notice that, although 
the gross Sip perturbation of Neptune grew much faster 
than Sip for Uranus, due to its higher k, the net pertur- 
bations for Neptune is considerably smaller, due to its 
low eccentricity. 

Later, in Section [Vl we will make a thorough best fit 
of these Sip perturbations to Keplerian orbits in order to 
find the minimal net perturbations for the three outer 
planets. 
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Uranus 
Neptune 




the plane of motion coincides with the plane given by 
6 = 7r/2. Writing the invariant g^uP^p" — — M 2 c 2 in the 
coordinates of eq. [15] and applying the constants of mo- 
tion of eq. [16] yields the following equations describing 
geodesic motion in a static, spherically symmetric gravi- 
tational field El: 



1 ( dr 
3 * dr 



dip _ I 

dr r 2 
dt _ £ 
dr~ B' 



(17) 
(18) 
(19) 



PIG. 4: The net perturbation to the argument of latitude 
after removal of the zero-order secular term 



III. PERTURBATIONS TO PLANETARY 
ORBITS - GENERAL METRIC THEORY OF 
GRAVITATION 

Now, let us generalize the calculation of the preceeding 
section and just assume that the gravitational field can 
be described by a metric theory of gravitation. This is a 
very general case that includes any PPN-parametrizable 
theory, including General Relativity. The purpose of this 
generalization is to include cases where the hypothetical 
perturbation to the gravitational field is not due to the 
presence of a gravitational source, but rather is caused 
by a deviation in the long-range behavior of the gravi- 
tational field from the behavior expected from General 
Relativity. As mentioned in section HI alternative gravi- 
tational theories have been proposed as explanations for 
the Pioneer anomaly, so we would like our analysis to be 
able to cover such cases. 



Here we have introduced the constants £ = Ej\ic 2 
(total energy per unit rest energy) and I = L/fi (to- 
tal angular momentum per unit rest mass). Notice that 
for bodies moving non-relativistically at speed v, £ — 

<Jb(1 + (v/c) 2 ) ~ 1 -I- i(f ) 2 - which is very close 
to unity. For the three outer planets, 1 — £ = 0(1O~ 10 ). 

Notice that these equations are valid for any metric 
theory of gravitation, and that we are not making any 
further assumptions on the gravitational interaction be- 
sides minimal gravitational coupling, i.e. a body moving 
freely in a gravitational field moves geodesically. 

In Newtonian gravity, a constant radial acceleration 
corresponds to an additional term in the gravitational 
potential that is linear in r. We then write the met- 
ric functions A and B as sums of the metric function 
B = 1 — 2GM/c 2 r of General Relativity (represent- 
ing the unperturbed, expected value of the gravitational 
field) and a small perturbation, representing the metric 
disturbance necessary to generate the Pioneer anomaly: 



A(r) = B (r) + a(r),B{r) = B (r) + f3(r). 



(20) 



A. Deriving the disturbance to the gravitational 
field 

Our first task is to analyze what constraints the in- 
troduction of a constant anomalous radial acceleration 
would put on the metric of the surrounding space. We 
assume a spherically symmetric metric, and begin by 
studying the dynamics of a body moving geodesically in a 
static, spherically symmetric, helio-centric gravitational 
field that can be described by a metric of the form 



ds^ 



-S(r)c 2 dt 2 +A(r) _1 dr 2 - 



2 (d(9 2 +sin 2 6>d^ 2 ) (15) 

An object with rest-mass \x > moving geodesically in 
this gravitational field has two constants of motion; the 
total energy, E, and the angular momentum, L: 



E 



, dt 
dr 



2 dip 
^ dr"' 



(16) 



where r is the proper time of the moving body. The 
motion is planar, and we have chosen coordinates so that 



The metric ds 2 , = -B c 2 dt 2 + B^dr 2 



r 2 df! 2 



sponds to a spherically symmetric vacuum gravitational 
field around a central body with mass M. We assume 
a and /? both to be much smaller than Bq, which allows 
us to expand perturbatively around well-known solutions 
for the spherically symmetric vacuum field case. We will 
refer to the metric ds 2 , as the unperturbed metric. Ex- 
pressing the radial equation, eq. [17\ to first order in a 
and j3, we get 



1 



drY 
dr" ) 



£^_ 

B 



1 a 



£■- 



Bo 
B a (l 







(21) 



In the following, we restrict our analysis to non- 
relativistic motion and weak gravitational fields. In this 
regime, the a-dependence of eq. [21] can be discarded, 
because, with v being the velocity of the moving object, 
■g 1 = (!) 2 <<1. Differentiating eq. [21] with respect 
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to proper time r, we get the expression for the radial 
acceleration: 



where prime denotes differentiation wrt. the radial co- 
ordinate r. The first term on the right-hand side of this 
equation is the expected radial acceleration generated by 
the unperturbed vacuum metric. The last term on the 
right-hand side can be seen as an anomalous acceleration 
produced by an anomalous term in the goo component 
of the metric. In order for — \c 2 j3' to be a constant, neg- 
ative acceleration, f3 must be linear in r: 



= Pr, 



(22) 



with P being a positive constant of dimension 1/r. It 
should be noted that in order to compute the anomalous 
acceleration observed by a distant observer, we should 
transform eq. [21] to the proper reference frame of the ob- 
server, which for a distant observer would correspond to 
coordinate time t at the location of the observer. How- 
ever, as |dt/dr — 1| = O(10 -9 )for the velocities and as- 
tronomical distances we are considering, we can safely 
equate proper time r with coordinate time t. If ap is the 
anomalous acceleration, directed towards the Sun, eq. [22l 
immediately gives us the relationship between the metric 
parameter P and ap: P = 2ap/c 2 £ 2 . 

We can conclude that an observed constant anomalous 
acceleration uniquely constrains the dominant part of the 
disturbance to the gravitational potential. An additional 
term, linear in r, must be present in the goo component 
of the metric in order to produce such an effect. This 
result is the same we get from Newtonian gravity. The 
analysis is independent of how this metric disturbance is 
produced, and is therefore also independent of the par- 
ticular gravitational theory applied to analyze the the 
disturbance. 



The unperturbed (a p = 0) equation has the following 
approximate solution valid for low eccentricity orbits [25. 



1 

uo 



ro 



a(l 



1 + ecos(/ (ip - cjq)) ' 



(25) 



where e is the eccentricity of the orbit, ujq is the ar- 
gument of the perihelion at a time when the planet is at 
perihelion, and f 2 — 1 — §{^f) 1 as r oo. / is 
very close to unity. When setting / = 1 we get the closed 
Keplerian orbits of eq. [4] In general, a factor / < 1 yields 
a small relativistic perihelion precession, which for low 
eccentricity orbits is ~ 6wGM / '(c 2 a(l — e 2 )) per cycle. 
Evaluating the relativistic perihelion precession for the 
three outer planets gives values of of the order of 1CT 3 
arc second per cycle. This is well below what would be 
observable, given the present data and their accuracy, 
which is of the order of 0.5 arc second. 

Let us proceed as we did for for the Newtonian case of 
Section ITlBl and assume that uo is a reference solution to 
the unperturbed equations of motion, of the form given 
by eq. [25] Let us then define the perturbation Su — 
u — uq, where u is a solution to the perturbed equation 
of motion, eq. EH having the same angular momentum 
i as uq. We obtain the following equation of motion for 
the perturbation 5u: 



d 2 5i 



(f 2 + 2k) 8u 



a(l - e 2 )(l + ecos(/^)) 



p (26) 



where v 



(f — LOo-, K 



"(I 



GM 



as m 



1 — f 2 is of the order 
Since K is e>(l(T 4 ), 



the Newtonian case before. Now, 
10 -9 for the three outer planets, 
we can safely ignore any Post-Newtonian corrections to 
the solutions, set / = 1 and proceed with the Newto- 
nian equation of motion, eq. [7] and the corresponding 
solutions, as derived in section ITT1 



B. Equations of motion 

Having found the perturbation to the gravitational 
field, we are now ready to explore the equations of mo- 
tion for a body moving geodesically in this field. When 
discarding the a-term of eq. [21] and using that j^f3 ~ p, 
we get the radial equation of motion 

By using eq. HH1 expressing the equation in terms 
of the reciprocal radial coordinate u = 1/r, taking the 
derivative of the equation with respect to ip and using eq. 
[22] we get the following second-order equation for u(ip): 

d 2 u ( 3GM \ GM a P 1 



IV. ANALYSIS METHOD 

A. Observations and their residuals 

The procedures that are applied for calculating the 
ephemerides of the planets use available theory, taking 
into account known gravitational sources, and carefully 
fitting the computations to available observational data 
[IB, E3|- One would therefore expect that any distur- 
bance to the gravitational field within the solar system 
would, if large enough, over time generate residuals in 
the position measurements of the planets as they would 
tend to drift from their expected positions. The question 
then is, what residuals would arise if we were to explain 
the Pioneer anomaly as the result of a perturbation to 
the gravitational field in the outer solar system, and how 
do they compare with the observed residuals? In this 
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TABLE I: Standard deviations in rii 


$ht ascension residuals for 


the three outer planets 




Planet 


Std. dev - ai (arc sec.) 


Uranus 


0.3 


Neptune 


0.3 


Pluto 


0.8 



section, we will elaborate the analysis method which we 
will apply in order to answer these questions. 

Before we start our analysis, let us review what the 
observables are and how to compare models with obser- 
vations. See Appendix [B] for a review of terminology and 
notation used. When restricted to performing optical ob- 
servations only of a distant object, which has been the 
case for the three outer planets, there are just two observ- 
ables: the right ascension (r.a.), and declination (dec). 
These are the position angles that locate the object on 
the celestial sphere. Accurate optical measurements of 
the positions of the outer planets have been performed 
since the early 20th century. The a priori uncertainty of 
these measurements were 0.5 — 1 arc second during most 
of this period In recent years, the a priori accuracy 
has been improved to 0.2 arc second. Computational 
models of the motion of the known objects in the Solar 
System can then be applied to compute the ephemerides, 
which are the expected orbits, of these objects. The free 
parameters of the model are fixed by making a best fit 
of the computations with the observations. The results 
are the Solar System Ephemerides, composed of highly 
accurate maps of the motions of the planets as well as 
other known Solar System bodies fl6l.ll7j|. 

Figure [5] shows the DE414 residualsflgj in right ascen- 
sion for the three outer planets. These residuals are the 
difference between the measured and ephemeris right as- 
cension values, taken over the entire span of available 
data. As could be expected, the residuals are smaller to- 
wards the end of the data span, apparently due to the in- 
creased accuracy of the measurement techniques applied. 
TableUgives the standard deviations of the r.a. residuals 
for the three outer planets[32|. The standard deviations 
of these residuals are within the a priori measurement 
accuracies that have been reported fl6Lll7j|. 
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FIG. 5: Actual r.a. residuals (in arc seconds) for the three 
outer planets, courtesy of E.M. Standish, JPL 



B. Description of the analysis method 

Our best option for determining whether or not a per- 
turbation to the gravitational field would have an ob- 
servable effect on the motions of the outer planets is to 
solve the basic equations of motion for the three planets 
under the hypothesis of a perturbed gravitational field. 
The purpose of these solutions, referred to as the simu- 
lated observations orbits, is to provide targets for best fit 
approximations using solutions to the unperturbed equa- 
tions of motion. Just like real ephemerides emerge as best 



fits of solutions of the equations of motion to the real ob- 
servations, we will construct simulated ephemerides by 
fitting solutions of the unperturbed equations of motion 
to the simulated observations orbits[33|. The simulated 
ephemerides obtained this way therefore play the same 
role as the real ephemerides, i.e. the best fit of known 
models to observations. From the simulated observations 
orbit and the corresponding simulated ephemeris, we can 
derive the corresponding values of the optical observ- 
ables; the right ascension a and the declination 8 [HiJ and 
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their residuals. These residuals (simulated obs. value - 
simulated ephemeris value), which are derived from the 
basic equations of motion, provide a realistic lowest or- 
der estimate of the residuals in the observables to be 
expected from a perturbation to the gravitatational field 
of the kind required to account for the Pioneer anomaly. 

Although we work with models that are idealized com- 
pared to the realistic models underlying modern plane- 
tary ephemerides, we still expect the simulated residual 
estimates to be realistic to lowest order in the pertur- 
bations, because to lowest order, perturbations decouple 
and can be studied independently. 

Finally, since we would like to compare the results of 
the hypothetical model to real observations, we will de- 
fine the predicted observations to be the sum of the real 
ephemeris values and the best fit simulated residuals. 
The underlying assumption here is that the simulated 
residuals can be used to predict the offset between the 
ephemeris and the expected observations under the hy- 
pothesis of a perturbed gravitational field. This is the 
closest we can get with respect to checking this hypoth- 
esis against data. We are not interested in the predicted 
observations themselves, only in their residual relative to 
the real observations. As we will see later, the (real - 
predicted) residual can be estimated by combining the 
real and simulated residuals. 

Having obtained the simulated residuals as the differ- 
ence between the simulated observations and the simu- 
lated ephemerides, the remaining problem is how to com- 
pare them with the real residuals, described in section 
IIV Al above. Let us formulate two hypotheses that we 
would like to test by statistical methods: 

HI: Given the uncertainty of the observations, the 
simulated r.a. observations do not deviate significantly 
from the simulated ephemerides 

H2: Given the uncertainty of the observations, the 
predicted r.a. observations do not deviate significantly 
from the real observations 

Hypothesis HI is testing, under the assumption of a 
perturbed gravitational field, whether the ephemerides 
still would provide good fits to the observations. If not, 
it would contradict the present state of matters, because 
the ephmerides do provide good fits to the data. 

Hypothesis H2 tests, again under the assumption of a 
perturbed gravitational field, whether the model predic- 
tions are consistent with the actual observations. 

Now, let us restate the original hypothesis stated in 
Section II Bl above: 

HO: The Pioneers move geodesically in a perturbed 
gravitational field 

The logical connection between hypothesis HO and the 
other two should be clear. HO is falsified if any of hy- 
potheses HI and H2 are falsified, because falsification of 
HI or H2 would imply that a perturbation to the gravita- 
tional field required to induce the Pioneer anomaly is not 
present. On the other hand, HI and H2 can obviously 
not be used to verify HO. 

Notice that the results will critically depend on judi- 



cial estimates of the measurement uncertainties. For a 
given planet I, we will use the standard deviation 07 of 
the r.a. residuals shown in Figure [5] as the best estimate 
of the uncertainty in the observations. If the ephemeris 
provides a good match to the observations, it is reason- 
able to assume that 07 will be close to the measurement 
uncertainty. If not, 07 would be higher than the mea- 
surement uncertainty. In either case, using o\ as a best 
estimate of the uncertainty in the measurements seems a 
reasonable, if not conservative, choice. 

In order to test hypotheses HI and H2, we will perform 
a chi square analysis with a confidence limit of 99%. Let 
us elaborate how we will apply this method in our case. 
Given Ni measurements (&i,j,tj) of the right ascension 
for planet I at measurement times tj, the real residuals 
Pl,j(tj) for planet I is the difference between the observa- 
tion ai : j and the corresponding ephemeris value ai(tj): 
pi,j(tj) = &i t j — ai(tj). These are the residuals that are 
plotted in Figure [5] above, which we will refer to as the 
real residuals. 

Given a solution to the equations of motion for planet 
I under the hypothesis that the gravitational field is per- 
turbed by a constant radial acceleration a p directed to- 
wards the Sun, we can derive the set of simulated right 
ascension values af (tj). We will refer to these as our 
set of simulated observations. The simulated ephemeris 
is the solution of the unperturbed equations of motion 
that provides the best fit to the simulated observations 
orbits, and from it we can derive a set of simulated right 
ascension ephemeris values af (tj ) . The simulated resid- 
uals pfj are then the difference between the simulated 
right ascension values and the corresponding simulated 
ephemeris values: pfj(tj) = &u(tj) ~ a f(tj)- 

The predicted r.a. value afj(tj) is the sum of the real 
ephemeris value ai(tj) and the simulated residual pf,: 

= aifa) + pfj (27) 

Define the real-predicted observation residual pf - to 
be the difference between the real observation and the 
predicted observation: 

pf'j = \j (tj) - "Cife) = Pufa) - Plj(*j) (28) 

We see from eq. [28] that we can estimate the real- 
predicted observation residuals by combining the real 
residuals pi,j(tj) with the simulated residuals pfj(tj). 
This allows us to compare the predicted observations 
with the real observations using the chi square method. 

Let af be the variance of the real residuals for planet I. 
As discussed previously, we will assume that 07 is a good 
estimate of the measurement uncertainty. Now, define 
the following chi square values that we will use to test 
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the hypotheses above: 

1 Nl 

^,i = ^E^(*i) a (29) 
i=i 

1 Nl 2 

&4^(^) (3o) 

i j=i 

Provided measurement errors are distributed normally, 
a model deviates from observations with probability p if 
the x 2 °f the residuals between the observations and the 
model prediction exceeds N[ + C p , where C p is the confi- 
dence limit applied. C p depends on the probability given 
(the confidence level) and the number of free parame- 
ters. In our case, the number of free parameters is 4, 
and we have chosen a confidence level of 0.99. This gives 
C p = 13. 

xfi can be used to test hypothesis HI above; if Xu > 
Ni + C p , the deviation of the simulated r.a. observations 
from the simulated ephemeris is statistically significant at 
a significance level of p. Similarly, \i i can be used to test 
hypothesis H2 above; if X12 > Ni + C p , the deviation of 
the projected r.a. observations from the real observations 
is statistically significant at a significance level of p. 



to the true anomaly, y, by tp = u> + v, where ui is the ar- 
gument of periapsis [33|. In our case, w is the argument 
of periapsis of the reference orbit, i.e. a fixed position. If 
ipa(t) is the argument of latitude of the reference orbit, 
ip p (t) = <fio(t) + Sip(t) is the argument of latitude of the 
perturbed orbit. Thus, we can derive expected values for 
the a and 5 observables from ip p by performing a coor- 
dinate transformation of the position vector (r p ,ip p ) in 
the orbital plane of the planet to geocentric equatorial 
coordinates. 

The second perturbation, which we will refer to as our 
simulated ephemeris, arises when we attempt to find 
the solution that best fits the simulated observations or- 
bit (r p , ip p ). In order to simplify the analytical treat- 
ment, we assume that this solution can be treated as 
a perturbation, denoted A<p, to the reference orbit as 
well. Having found the best fit solution by performing 
a least squares approximation of the second perturbed 
orbit (r E (t) = r (t) + Ar{t),cp E (t) = <po(t) + A<p(t)) to 
the simulated observations orbit (r p , cp p ), we can estimate 
the simulated residuals in a and 5 as the difference be- 
tween the values of the simulated observations and the 
simulated ephemeris, as described in section ITVBl 



V. RESIDUAL ANALYSIS 



C. Comments on the perturbative approach for 
solving the equations of motion 

In this paper, we study what observational effects a dis- 
turbance to the gravitational field would have on orbits of 
a planet moving in a static, spherically symmetric grav- 
itational field around the Sun. We work perturbatively 
to lowest order in the perturbation. Since we are looking 
for dominant effects only, we are free to disregard other 
perturbing effects, such as the presence of other bodies in 
the Solar System. These simplifications allow us to treat 
the problem analytically. 

In section [Ul we started our analysis by choosing a 
reference orbit, which is an idealized orbit satisfying the 
unperturbed equations of motion and matching the or- 
bital parameters of one of the three outer planets. In 
the Newtonian case, the reference orbit is a Keplerian 
orbit that is set up by applying the mean Keplerian ele- 
ments of the planet as initial conditions. This reference 
orbit has no other role than being the reference for two 
perturbations: 

In section |Ul we found the first pertubation, which 
gives rise to the simulated observations. It comes from 
solving the equations of motion for a perturbed gravi- 
tational field. For this orbit, the simulated observations 
orbit, we applied the assumption that the gravitational 
field is perturbed in a way that results in a constant, 
anomalous radial acceleration ap = 8.7 x lCT 10 m/s to- 
wards the Sun. By integrating the equations of motion, 
we obtained the perturbation to the argument of latitude, 
denoted by Sip. The argument of latitude, ip, is related 



If (ro(vo)) Vo(*)) is the reference orbit, (r p (ip) = 
ro(<Po) + Sr(<p ),ip p (t) = ipo (t) + 5ip(ipo)) is the simulated 
observations orbit, i.e. the perturbed orbit from which we 
can derive simulated observations. The simulated obser- 
vations orbit was derived in section [TTJ Now it is time to 
proceed with the next step of our analysis: Following the 
analysis procedure outlined in Section IIVBI above, our 
next task is to compute the simulated ephemeris, which 
is the Keplerian orbit that provides the best fit to the 
simulated observations orbit. 



A. Finding the simulated ephemerides 

Eq. [3] represents the unperturbed equation of motion 
that we will solve. Since the simulated observations orbit 
(r p (ip), <p p (t)) is just a perturbation of the reference orbit 
(ro, ipo), we can assume that the best fit Keplerian orbit, 
our simulated ephemeris, also represents a pertubation to 
the reference orbit. Perturbing eq. [3] gives the equation 
for the radial perturbation Au = u — uq, where u here 
is another Keplerian orbit with angular momentum t = 
In + A£: 



d 2 Ai 



Au = -- 



2 A£ 



d<p z o(l - e 2 ) e 

Eq. [31] has the general solution 



(31) 



If 2 A A 

Au=- eCisin v + eC 2 cos v - — oK~r~ ] , (32) 

a V I 1 _ e ) *o / 
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where C\, Ca and Ai will be treated as free parameters. The exact solution to eq. [33] is 

Similarly to how we derived eq. HH we get the equation 
for the angular perturbation: 

dAw Ai Ar , , 

d</5 £ T-o 

I 



Av? = Co - 



A£ 



2C 2 (1 



2Ci(l-e 2 )ln(l + ecos v)- 



s/T 



Ai 



arctan( 



e v . 
- tan 2 } 



r 



J341 



Define uk = uo + Au and i^k = <^?o + Ay>. Ay; in- 
troduces another free parameter Co- This leaves us with 
four free parameters that we will fit to the simulated or- 
bit: C , Ci, C 2 and At. 

Define the simulated residual Aip by Aip = ip p (ip ) — 
VkCvo) = Sip — Aip. The simulated residual Aip can be 
written in the form 

Ai 

Aip = Sip - C + T fo(v>) + Ci/ifa) + Ca/afoO, (35) 

where /o, /1 and / 2 are functions that to third order in 
eccentricity take the form 

fo = (3 + 2e 2 ) 1/ - e (4 + 3e 2 ) sin 1/ + e 2 sin 2v 

e 3 

sin3z/ + C(e 4 ) 



cos 2v 



6 

+ C(e 4 ) 



cos 3^ 



residual Aip can be eliminated by adjusting the angular 
momentum, the Pioneer effect is in principle still observ- 
able from the perturbation to the observation angle re- 
sulting from the non- vanishing radial residual Ar. In the 
following, we will make a best fit of the simulated angular 
residual Aip and ignore the effect of the remaining radial 
perturbation on the observable angle. 

Now, let us define the best fit solution (te, <Pe) as the 
solution {rx, <-Pk) that minimizes the variance of the sim- 
ulated residual Aip over the measurement interval [to, £1]: 



1 



v(ti) 



(ip(h) - ip(t )) 



dip Aip 



(36) 



¥>(*o) 



The best fit is made by finding the values of the free 
parameters At, Co , C\ and C2 that give minimal variance 
tr 2 . Equation [36] can be integrated analytically, and gives 
a second order polynomial in the free parameters Co, 
Ci, C2 and Ai. This polynomial has a unique minimum 
that can be found by solving the equations -St- c 2 = , 



= 0. 



a 

i)C 2 



ai = and 



0Ae u p 



= 0. 



ev - 



e . e . _ 

— smzi/ smo^ 

2 6 



+ G(e 4 ) 

Let us consider for a moment the circular limits of eqs. 
[321 and [35] In that case, fx = f 2 = 0. Eqs. M and [32] 
show that in the circular limit, the simulated radial resid- 
ual 8u — Au can be eliminated by properly adjusting the 
angular momentum of the the Keplerian orbit. However, 
in this case we are left with a non- vanishing simulated an- 
gular residual Aip. If we, on the other hand, adjust the 
angular momentum and Co to make the angular residual 
Aip vanish, we are left with a non-zero radial residual 
Ar = Sr — Ar = na/3. Thus, even in the circular limit, 
we can not simultaneously fit both the angular and ra- 
dial perturbations by a Keplerian orbit. When projecting 
Ar onto the celestial sphere for the near circular orbits 
of Uranus and Neptune, it translates to perturbations 
in the observation angles of 0.2" and 0.3", respectively. 
This implies that in the circular limit, even if the angular 



B. Analyzing the simulated residuals 

Figure [6] plots — pP , the predicted-real r.a. residuals, 
and the simulated r.a. residuals for the three outer plan- 
ets. In order to make the simulated residuals easy to 
compare with the real residuals, a coordinate transfor- 
mation of the orbit to the mean equatorial plane has 
been done. The plots also show the standard deviations 
of the real residuals. 

From Figure [6] we can immediately see that the simu- 
lated residuals for Neptune are very small compared to 
the real residuals. This implies that a Pioneer- type per- 
turbation would have no observable effect on the orbit of 
Neptune. Furthermore, by visual inspection of the real- 
predicted residuals in Figure [6] we see that there are 
significant deviations between real and predicted obser- 
vations for Uranus and Pluto. 

Now, let us put this on formal grounds and test hy- 
potheses HI and H2 formulated in Section IIV Bl above. 
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FIG. 6: Simulated (solid lines) and predicted-real r.a. residu- 
als (in arc seconds) for the three outer planets. Dashed lines 
indicate the standard deviation of the real residuals. 



Given the simulated residuals p p and the real residuals 
p, we must compute the chi square values x? i an d X12 
that target these hypotheses. Table ITTI shows the results 
of this computation and the tests of hypotheses HI and 
H2 for each of the three outer planets. 

The main conclusions drawn from this analysis are 
that, given an uncertainty ai in the r.a. observations: 

1. For Uranus, the simulated r.a. observations show 
a deviation from the simulated ephemeris that is 



TABLE II: Chi square values of the simulated residuals for the 
three outer planets and the results of the tests of hypotheses 
HI and H2 of Section HVBl 



Planet 


<*(") 


N t 


xh - N t 


xh ~ N i 


HI 


H2 


Uranus 


0.283 


3678 


486 


4930 


False 


False 


Neptune 


0.293 


3800 


-3799 


6.5 


True 


True 


Pluto 


0.771 


2119 


-1510 


773 


True 


False 



statististically significant well beyond the 99% con- 
fidence level. For Pluto, however, no such deviation 
is found for the simulated r.a. observations. Thus, 
HI is falsified for Uranus and verified for the other 
two planets. 

For both Uranus and Pluto, the predicted r.a. ob- 
servations show a deviation from the real r.a. ob- 
servations that is statististically significant well be- 
yond the 99% confidence level. Thus, H2 is falsified 
for Uranus and Pluto, but verified for Neptune. 



VI. CONCLUSIONS 

Under the assumption that the observed frequency 
drift of the Pioneer spacecraft is caused by a constant, 
anomalous acceleration that is directed towards the Sun, 
and that the Pioneer spacecraft move geodesically in a 
perturbed gravitational field that is static and spheri- 
cally symmetric, we have analyzed the effects that such a 
perturbation to the gravitational field would have on the 
orbits of the three outer planets. The basic assumptions 
are that the gravitational field can be modeled within a 
metric theory of gravity and that particles couple mini- 
mally to the gravitational field, i.e. they move geodesi- 
cally when falling freely in a gravitational field. Our re- 
sults are model independent within this class of theories 
in the sense that they do not depend on the particu- 
lars of the gravitational theory, e.g. the nature of any 
gravitational source or how gravitational sources couple 
to the spacetime metric. The analysis is valid in the 
weak field limit of any metric gravitational theory, and 
it therefore applies to any gravitational theory that can 
be parametrized within the PPN framework. 

We showed in section IIII AI that the metric perturba- 
tion that would be needed to induce the Pioneer anomaly 
can be uniquely constrained to lowest order. This im- 
mediately implies that perturbed orbits due to the per- 
turbation in the gravitational field will be the same, to 
lowest order in the perturbation, for all metric theories 
of gravity. Furthermore, by solving the basic equations 
of motion for a perturbed gravitational field, we indeed 
found that Post-Newtonian effects were small compared 
to Newtonian effects. 

By solving the unperturbed equations of motion and 
finding the solutions that provide the best fit to the per- 
turbed solutions, we were able to derive simulated right 
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ascension residuals that could be compared with the real 
ephemeris residuals. 

We then formulated an analysis method, based on the 
chi square method, for comparing simulated residuals 
with the real ephemeris residuals. We formulated two 
hypotheses that, if falsified, would falsify the hypoth- 
esis of the presence of a perturbation to the gravita- 
tional field. Following from Einstein's equivalence prin- 
ciple, this would imply that the Pioneers do not move 
geodesically in a perturbed gravitational field. Hypothe- 
sis HI tests, under the assumption of a perturbed grav- 
itational field, whether the ephemeris still provides a 
good fit to the observations expected from a perturba- 
tion to the gravitational field. If not, it would contradict 
the present state of matters, because the ephemerides 
do provide good fits to the data. Hypothesis H2 tests, 
again under the assumption of a perturbed gravitational 
field, whether the model predictions emerging from the 
model of a perturbed gravitational field are consistent 
with the actual observations. Falsification of either HI 
or H2 would imply that a perturbation to the gravita- 
tional field required to induce the Pioneer anomaly is 
not present, which under Einstein's equivalence principle 
would imply that the Pioneers do not move geodesically. 

The analysis shows that the ephemerides provide good 
fits to the perturbed orbits expected from a perturbation 
to the gravitational field for Neptune and Pluto, but not 
for Uranus. Thus hypothesis HI is falsified for Uranus. 
The fact that no significant deviation in the simulated 
residuals for Pluto is found can be attributed to the fact 
that the measurement uncertainty er; for Pluto is higher 
than for the other two planets. Furthermore, the analy- 
sis also shows that the deviation between real and pre- 
dicted r.a. observations is statistically significant for both 
Uranus and Pluto, implying that hypothesis H2 is falsi- 
fied for both planets. The statistical significance of this 
result is well beyond the 99% confidence limit in both 
cases. For Neptune, on the other hand, we find no signif- 
icant deviation, so hypothesis H2 is verified for Neptune. 
This is not surprising, because Neptune's orbit has very 
low eccentricity, and, as we saw from the discussion of 
this subject in Section IV A| the angular residual arising 
from the Pioneer effect can be fitted out in the circu- 
lar limit (the small perturbation to the observation angle 
that comes from projecting the remaining radial residual 
onto the celestial sphere is ignored in this study). 

It should be stressed again that these conclusions rest 
on judicial esitmates of the uncertainties in the observa- 
tions. We have chosen to use the standard deviations 01 
of the real residuals as the best estimates of the measure- 
ment uncertainties. 

Based on the available observations of the three 
outer planets, an assessment of their uncertainties 
and the available residuals between observations and 
ephemerides, we are therefore lead to the conclusion 
that the presence of a perturbation to the gravitational 
field necessary to induce the Pioneer anomaly is in con- 
flict with available data. This implies that the Pioneer 



anomaly can not have gravitational origin and, conse- 
quently, the Pioneer spacecraft do not follow geodesic 
trajectories through space. Hence, any external interac- 
tion model introduced to explain the Pioneer anomaly 
should, in one way or another, violate the equivalence 
principle. It should be stressed that "gravitational ori- 
gin" in this context really means any physical interaction 
that is mediated indirectly via the space-time metric only. 
It is assumed here that particles move geodesically when 
moving under gravitational influence only (minimal cou- 
pling), so models implying other, non-minmal coupling 
schemes can not be excluded by our results. Our conclu- 
sion is consistent with recent findings of Iorio and Giudice 
[20| as well as Izzo and Rathke [li}|. It should be noted 
though, that while Iorio and Giudice find large secular 
effects on the orbits of the three outer planets, we find 
no such effects. The reason for this discrepancy is that, 
as shown here, any secular effects would be canceled by 
the simulated ephemeris, i.e. the best fit solution to the 
unperturbed equations of motion. The only observable 
effect is therefore the net residuals that still remain after 
the ephemeris cancels the big secular deviations. 

This result can be used to rule out any model of the 
Pioneer anomaly that implies that the Pioneer space- 
craft move geodesically through space and that explains 
the Pioneer anomaly as the effect of a disturbance to 
the spacetime metric, regardless of how this disturbance 
is created. This includes some of the proposed models 
of the anomaly mentioned in section HI such as scalar 
field models [Icj], braneworld scenarios [ll[ and mod- 
els involving alternative metric theories of gravitation 
0, M, Ell, Q such as TeVeS [H] • 

We conclude that, in order to resolve the enigma of the 
Pioneer anomaly, there is a need for theories and models 
that may explain the apparently non-geodesic motion of 
the Pioneer spacecraft. 
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APPENDIX A: PERTUBATIVE EXPRESSIONS 

Here, we provide the expansions to third order in ec- 
centricity of the perturbative expressions given in Section 
III Bl above. All pertubative calculations were made with 
these third order expressions, so the results given in this 
paper are indeed valid to third order in eccentricity. 



15 



The C(e 3 ) perturbation in reciprocal radial distance, 
8u, which to lowest order is given by eq. is 



k f 20 o ^ e 
ou = — 1 + — e — — cos 2v + — cos "iv 
a V 8 2 8 

e A 5 2 \ 
1 + -e cosj/ 

«\ 2) 

+ A sin(Vl + 2nv) + B cos(Vl + 2k i/) + C(e 4 ) (Al) 

Demanding <5u = for k = fixes A and -B to take the 
values A = and B = f (1 + |e 2 ) + C(e 4 ). This gives an 
angular perturbation 



k f 20 2 e 2 „ e 3 

ou = — < 1 H e cos 2v -\ cos 6v 

a \ 8 2 8 

-e(l + ^e 2 )zysin^j +C(e 4 ). (A2) 



account of this terminology, see any astrodynamics text- 
book, such as [29j. 

TABLE III: Data used in the calculations 



Variable 


Uranus 


Neptune 


Pluto 


to 


1914-07-08 


1913-12-28 


1914-01-23 




06:59:46 


06:41:17 


18:58:11 


h 


2006-09-30 


2006-09-30 


2006-08-26 




05:45:39 


04:11:14 


02:45:13 


^o(deg) 


-220.57 


72.49 


-133.75 


fi(deg) 


171.85 


273.16 


40.20 


.Mo (deg) 


-224.19 


71.53 


-109.47 


.Mi (deg) 


171.06 


274.17 


24.74 


e 


0.0473 


0.00882 


0.251 


a (AU) 


19.19 


30.07 


39.53 



The third order expression for the radial perturbation 
Sr, corresponding to eq. [101 is 

Sr = —Ka jl + 2e 2 -e ^2 + ^e 2 ^j cosv + e 2 cos2zy 

3 5 
— -e 3 cos Sv — e(l + -e 2 )u sin v 
8 4 

+e 2 v sin 2v - ^e 3 ^ sin + C(e 4 ). (A3) 

Expanding the equation for the angular perturbation 
to third order in eccentricity gives 

^ = 2K-h+2e 2 -e(l + 2e 2 )cosi/ 
dip L 

e 2 1 
—ev{\ + 2e 2 — e cos u + — cos 2v) sini/ + -e 3 cos 3z/ 
2 8 

+ G(e 4 ). (A4) 

Integration of eq. IA4I wrt. ip gives the angular per- 
turbation corresponding to eq. HU 



Stp = K |i/(2 (l + 2e 2 ) + | (4 + 7e 2 ) 



(3 cos 2v — ecos3^^ 



6 6 

■- sin;/ (8 + 15e 2 - ecosi/) H sin 3^ ^ + C(e 4 ). 

36 1 



(A5) 



APPENDIX B: TERMINOLOGY AND 
NOTATION 

Here, we provide a short review of the astonomical ter- 
minology and notation used in the paper. For a complete 



The Ecliptic plane is the plane in which the Earth 
moves around the Sun. The Ecliptic is the cross section of 
the Ecliptic plane with the celestial sphere. The right as- 
cension (r.a.) and declination (dec.) angles fix the object 
on the celestial sphere. These angles are measured rela- 
tive to the celestial equator, which is the cross section of 
Earth's equatorial plane with the celestial sphere. The 
right ascension, denoted a, measures the angle of the 
object eastward along the celestial equator from a fixed 
point on the sky; the vernal equinox, denoted T, which 
is the point where the Sun crosses the celestial equator in 
the spring. The declination, denoted S, measures the an- 
gle northward from the celestial equator. The ascending 
node is the direction where the object ascends through 
the ecliptic. The argument of latitude (here denoted ip) 
is the angular distance of the object measured in the or- 
bital plane from the ascending node, along the direction 
of movement. The argument of periapsis (denoted w) is 
the angular position of the periapsis, which is the point 
in the orbit where the object is closest to the Sun. The 
true anomaly, denoted v, is the angular position of the 
object, measured from the periapsis, in the orbital plane 
and in the direction of movement, v = tp — lo. 



APPENDIX C: DATA USED IN CALCULATIONS 

Table IIIII lists various constants used in the residual 
calculations. 

to = Time of initial data, t\ = Time of final data, 
1/q = v(to) = true anomaly at time of initial data (deg), 
V\ = v{t\) — true anomaly at time of final data (deg), 
Mo = M(to) = mean anomaly at time of initial data 
(deg), Mi = M(ti) = mean anomaly at time of final 
data (deg). 
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